Self-organized criticality in linear interface depinning and sandpile models 
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The dynamics of an elastic interface profile h(x, t) under a driving force increasing at rate c, a 
restored force —eh, and disorder is investigated. Using perturbation theory and functional renor- 
malization group the phase diagram and the scaling exponents, up to the first order in e — 4 — d, are 
obtained. The model is found to be critical in the double limit e — » and c/e — > and belongs to 
a different universality class as that of constant force models. It is shown that undirected sandpile 
models with stochastic rules and linear interface models with extremal dynamics belong to this new 
universality class. 
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The problem of interface roughening in the present of 
quenched disorder is a topic of recent interest, due to its 
importance as a paradigm in condensed matter physics 
and due to the broad range of applications. In general a 
d-dimensional self-affine interface, described by a single- 
valued function h(x,t), evolves in a (d + l)-dimensional 
medium. Usually some form of disorder affects the mo- 
tion of the interface leading to its roughening. Earlier 
studies jlj focuses on time-independent uncorrelated dis- 
order but most recent studies analyses the motion of in- 
terfaces under quenched disorder In the present 
of quenched disorder an constant force driving two uni- 
versality classes has been found ||. One is described 
by the Kardar-Parisi-Zhang equation jl| with quenched 
noise. In this case the interface is pinned by paths on a 
directed percolation cluster of pinning sites || . The sec- 
ond class is described by the Edwards- Wilkinson (EW) 
equation with quenched noise, linear interface depin- 
ning (LID). 

The interface may be driven either by extremal dynam- 
ics — |l"o|] or by constant force j|Tl|— ^L3|| . While constant 
force models has been extensively studied in the litera- 
ture, either by functional renormalization group (FRG) 
jl4],[l5| or numerical simulations (lj] , extremal dynamics 
models are less known. However, in the last years ex- 
tremal dynamics models have gain more attention due 
to its relation with the theory of self-organized criticality 
(SOC). 

SOC was introduced to explain the critical behavior of 
a vast class of driven dissipative system which evolve into 
a critical state . In its early state it was believed that 
such a critical state is insensitive to changes in control 
parameters and no fine-tuning is needed. More recent 
interpretations of this phenomena have shown that criti- 
cality in SOC systems is obtained after some control pa- 
rameters, for instance the driving and dissipation rates, 
are fine-tuned to zero [JlTj . 

In the present work we focus our attention in the mo- 
tion of a d-dimensional interface profile obeying the fol- 
lowing equation, 



\dth = TV 2 h + F - eh + ri{x,h). (1) 

with a force increasing at constant rate, F — ct. Here A 
is a friction coefficient, T is the surface tension, c is the 
driving rate, and e a nonnegative constant. The random 
force r]{x, h) is Gaussian distributed with zero mean and 

(r)(x, h)n{x', ti)) = S d (x - x')A(h - ti), (2) 

where A(/i) is a monotonically decreasing function. 

The aim of this work is to show that models described 
by this equation are in a different universality class as 
constant force models, and that undirected sandpile mod- 
els with stochastic rules and LID models with extremal 
dynamics belongs to this new universality class. The mo- 
tivation for the analogy between sandpile models and LID 
models is based in a work by Paczusky and Boettcher 
|]i~8| , where it is shown that a one dimensional critical 
slope sandpile model is in the same universality class as 
the depinning transition of a d + 1 interface dragged at 
one end. Although this work pointed out the analogy 
between LID and sandpile models its analysis was lim- 
ited to one dimension. We onclude that the existence of 
"spontaneous" criticality in LID with extremal dynam- 
ics, as in sandpile model Jl7t , is just a consequence of the 
unprecise definition of these models. 

If the force F is constant and e = then eq. (|l|) 
is reduced to the EW equation with quenched noise. 
This case has been extensively studied in the literature 
|l4],[l5| . A depinning transition takes place at certain crit- 
ical field force F c determined by the disorder. For F < F c 
the interface is pinned after certain finite time while for 
F > F c it moves with finite average velocity which scales 
as v ~ (F - Fc)' 3 . 

When e > and the force increases at rate c then the 
interface is never pinned by disorder, but always moves 
with a finite average velocity v. A perturbative solution 
of eq. (|l|) can thus be found expanding h(x, t) around the 
flat co-moving interface vt. Taking h(x, t) — vt + y(x, t) 
we obtain the following equation for y(x,t) 

\d t y = TV 2 y + (c - ev)t - eh - \v + r](x, vt + h). (3) 
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The average velocity is obtained using the constraint 
(y(x, t)) — 0. For this purpose is better to work with the 
equation for the Fourier transform h(k,w). The effective 
external held (c — ev)t gives a singular term of the or- 
der of lo~ 2 . This singular term predominates in the low 
frequency limit resulting, after imposing (y(k,uj)) = 0, 



(4) 



which is valid to all orders of perturbation expansion. 

Another exact result can be obtained if one com- 
putes the low-frequency and long- wavelength susceptibil- 
ity. Adding a source term tp(x, t) to the right hand side 
of eq. (||) and going to the Fourier space one obtains the 
generalized response function 



G(k,w) 



where 



h(k, uj) 



<p(k,ui) 



=0 / [G (k,u;)}-i-Z(k,u)' 



[Goik^uj)}- 1 =Tk 2 -iXoj + e 



(5) 
(6) 



is the bare correlator and £(fc,w) is the "self-energy". 
Since £(0,0) = and Go^O)" 1 = e it results that the 
low- frequency and long- wavelength susceptibility (or sim- 
ply the susceptibility) is given by 



X = G(0,0) 
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This result is also exact to all orders of perturbation ex- 
pansion. Thus, when e — > the susceptivility diverges 
and, therefore, the system is critical. 

To go further we perform a FRG analysis of the prob- 
lem, following the general ideas developed for the con- 
stant force case ]l5| |. We construct the generating func- 
tional Z = J DhDhexp(S) with action 

S = J d d xdtih[\d t h - TV 2 h - F + eh - rj(x, h)], (8) 

where h(x, t) is an auxiliary held. After averaging over 
disorder an expansion around the mean-held (MF) solu- 
tion yields a generating functional with the low frequency 
form Z J DyDye~xp(S), where the effective action S 
is given by 



S = - I d d xdt{[F - F MF (v)]y(x, t) + y(x, t){Xd t - TV 2 + e)y(x, t)} 
d d xdt 1 dt 2 y(x, h)y(x, t 2 )C[v(t 1 - t 2 ) + y(x, ti) - y(x, t 2 )], 



(9) 



where y and y are coarse grained versions of h — vt and 
— ih, respectively, and C(h) is the MF correlation func- 
tion. Two differences appear with the constant force case. 
First, in our case F = ct and Fmf(v), the MF force cor- 
responding to a velocity v, is given by Fmf(v) — cuF(v)t. 
Moreover, we have obtained that c — ve (see eq. (^)) ex- 
actly, within and beyond the MF approximation. Hence 
F — Fmf(v) = 0. Second, in the Gaussian part of the ef- 
fective action there is an extra term associated with the 
restored force, characterized by the coefficient e. As it 
was shown above the susceptibility diverges when e — > 0. 
Thus, e is the control parameter of the interface described 
by eq. (0). On the contrary in the constant force case 
e = and F — Fmf(v) is the control parameter. 

The RG transformations are carried out as follows. We 
integrate out the degrees of freedom in a momentum shell 
near the cutoff A and rescale x — > bx, t — > b z t, y —> b^y, 
and y — > b e ~ d y, where b = e l with I — > 0. As usual the 
cutoff A appears because we start our analysis from a 
coarse-grained equation, where we cannot resolve spatial 
details smaller than A -1 . 

The renormalization of the e term yields 

d ± = (e + z + Qe, (10) 



which implies that the correlation length scale as £ ~ e v 
with v = 1/(9 + z + (). Since 9 + z + C) = 2 ff§ we fi- 
nally obtain v = 1/2, which differs from the one obtained 
in the constant force case, due to the existence of differ- 
ent control parameters. On the contrary other scaling 
exponents results identical. For instance |lj| 

c = |, .->-§.. (ID 

On the other hand, vt and y must scale in the same, so 
that v — > b^~ z v yielding 

f = (C"*>- (12) 

Thus, to reach the critical state both e and v should be 
hnc-tunncd to zero, i.e. criticality is obtained in the dou- 
ble limit e — > and v = c/e — > 0. From eq. ( |l2] ) we define 
the characteristic velocity v c ~ with 

(3 = v(z-Q. (13) 

Note that in constant force LID the average interface 
velocity in the supercritical regime is given by v ~ 
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(F — F C )P. On the contrary, in the present model the 
average interface velocity is fixed through eq. (||) and 
the exponent (3 just characterizes the e dependence of 
the characteristic velocity v c , which delimit the subcriti- 
cal and supercritical regimes. For v 3> v c the noise term 
in eq. (^)is approximately annealed obtaining the EW 
equation with annealed noise. In this case ( = for 
d > 2. The driving field predominates over disorder and, 
therefore, the model is supercritical. 

In the subcritical regime e > and v -C v c the dynam- 
ics takes place in the form of avalanches, characterized by 
the avalanche size distribution P(s) = s~ T g(s/s c ), where 
s c ~ e" 1 /"' is a cutoff avalanche size. In the critical state 
s c L D , where D is the avalanche dimension, and £ ~ L 
leading to the scaling relation a = 1/Dv. Another scal- 
ing relation is obtained taking into account that \ — ( s }i 
leading to 7 = (2 — t)/o. On the other hand, for d < d c , 
the avalanche dimension and the roughness exponent are 
related via D = d + ( JlIJ . Using these scaling relations 
and the values for 7, u, £, and z computed above we 
obtain 

D = d + (, r = 2(l-L>- 1 ). (14) 

To investigate the analogy with sandpile models let us 
analyze a discretized variant of eq. (Q). In the cellular 
automaton version of eq. (Q) one defines the total force 

Fi=J2 H J - 2dH > + ht - eH > + ^(^)> ( 15 ) 

nn 

and sites where Fi > are updated in parallel, advancing 
the interface by one Hi — > Hi + 1. Here ^2 nn Hj — 2diJ; is 
a discretized Laplacian, where the sum runs over nearest 
neighbors. Instead of follow the evolution of the interface 
profile Hi one may keep track of the total force Fi . The 
evolution rules for the total force Fi are given by: 1- on 
each step Fi is increased by h in all sites and 2- all sites 
where > are updated in parallel according to the 
toppling rule 

Fi -> Fi-2d-e + rjai , , 

F — > F + 1 (16) 

± nn T ± nn \ x > 

where r\ a i = rji(Hi + 1) — rji(Hi) is a zero mean uncorre- 
cted annealed noise and nn denotes nearest neighbors. 

From a simple inspection of eq. (^) one can see that 
the total force follows the evolution rules of an undirected 
sandpile automaton with an annealed noise, under a driv- 
ing field h and with dissipation rate per toppling event 
e. This class of sandpile models has been studied by 
Vespignani and Zappery Jl7[ , using mean-field and field 
theories. They have obtained that 7 = 1 and v — 1/2 are 
exact in all dimensions. Moreover, their field theory pre- 
dicts an upper critical dimension d c = 4. These results 
are in agreement with those obtained here and strongly 
suggest that both models are in the same universality 
class. 



In the mean-field and field theory by Vespignani and 
Zappery conservation (understanding conservation as the 
balance between input and output energy) is a necessary 
condition for stationarity, which implies that the density 
of toppling (active) sites is given by p a = h/e. The equiv- 
alent condition in our approach is found in eq. , which 
is a consequence of the balance between the driving force 
ht and the average restored force evt. The balance of 
these two forces is thus our necessary condition for sta- 
tionarity. Moreover, since the interface advances only 
at active sites (those where Fi > 0), in one unit, then 
v = p a , making the connection between our approach 
and that of Vespignani and Zapperi. Is also important to 
note that v < 1, which is consistent with the fact p a < 1. 

Now, we proceed to show that the extremal LID model 
corresponds to the critical state of eq. (|l|), i.e. h, e — > 
and v — ► 0. The condition h — * carry as a consequence 
that, if at time step to there are no active sites (F{ < 
for all sites) then at time step t\ = —Fj(to)/h the site 
j, with the maximum total force, will be active. In the 
language of interface depinning the interface at site j will 
advance one unit, in the sandpile language the site j will 
topple. It is thus clear that if no active site is present 
the system will follows extremal dynamics. Now, what 
happened during the evolution of an avalanche? 

The active site j will transfer energy to its nearest 
neighbors, which at the same time may become active, 
and so on, an avalanche is generated. It is thus possible 
that at certain time step t there will more that one active 
site. These sites will be updated in parallel according to 
the evolution rules described above. However, the order 
in which these sites are updated is not important. The 
process of toppling can never transform any active site, 
different from itself, in inactive and, therefore, the others 
active sites will remain active. On the other hand, the 
energy transferred to its neighbors is constant, indepen- 
dent of the total force at this site. Moreover, the site 
with the maximum total force will be always among the 
set of active sites. Hence, one can arrange the sequence 
of toppling events, of actives sites at time step t, in such 
a way that always the site with the maximum total force 
topples. 

Finally, in the double limit h — > 0, e — > the total 
force is reduced to the local force 

nn 

Hence, in the critical state of the cellular automaton ver- 
sion of eq. (^) one may arrange the sequence of toppling 
events in such a way that the site which topples has the 
maximum local force /j, which is the extremal dynamics 
variant of LID. We thus conclude that extremal dynam- 
ics models corresponds to the critical state of the LID 
model described by eq. (|l|). 

In table | some numerical estimates for the LID model 
with extremal dynamics and the Manna d-state model, 
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the prototype of stochastic sandpile model, are given. 
Our FRG estimates from eqs. ( |TT| ) and ([u]) are also 
shown for comparison. Increasing the lattice dimension 
the FRG predictions get closer to the numerical esti- 
mates, obtaining a complete agreement in three dimen- 
sions. Using the numerical estimates for r and D, re- 
ported for the LID model with extremal dynamics, we 
have tested the scaling relation (2 — t)D — 7Z/~ 1 = 2, 
within the numerical error it is fulfilled. In all cases the 
difference in the numerical estimates for both models are 
contained in the error bars, which implies that they are 
in the same universality class. Our predictions are thus 
confirmed with the numerical simulations. 

In the discussion we have not included the Bak-Tang- 
Wiesenfeld (BTW) model neither a sandpile model with 
stochastic dissipation introduced by Chessa et al [^9| . 
In the first case because the BTW model is determinis- 
tic and it is not clear yet if stochastic and deterministic 
models are in the same universality class. The LID model 
introduced here was mapped into a sandpile model with 
annealed noise and, therefore, our conclusions cannot be 
extended to deterministic models. In the second case be- 
cause according to the simulations by Chessa et al JToj ] 
the upper critical dimension for their stochastic model 
is 6, in disagreement with d c — 4 as obtained from our 
FRG analysis and the numerical simulations of the d- 
state Manna model ^O). On the theoretical side, d c = 4 
is in agreement with previous reports by Dfaz-Guilera, 
using a dynamic renormalization group approach for the 
BTW and Zhang models EL and by Vespignani et al 
& 

In summary, we have shown that the existence of 
"spontaneous" criticality in LID with extremal dynam- 
ics, as in sandpile model Jl7| , is just a consequence of the 
unprecise definition of these models. SOC corresponds 
to the onset of nonlocality in the dynamics of the inter- 
face. Nonlocality, and hence criticality, is obtained by 
fine tuning the control parameters, precisely as in con- 
tinuous phase transitions. The extremal dynamics corre- 
sponds with a fine tuned interface depinning transition 
at constant velocity. It was also demonstrated that LID 
with extremal dynamics and undirected sandpile mod- 
els with stochastic rules belong to the same universality 
class, which is different from that of constant force LID. 

We thanks A. Vespignani for helpful comments and 
suggestions. This work was partially supported by the 
Alma Mater prize, given by the University of Havana. 
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TABLE I. Scaling exponents for the LID model with ex- 
tremal dynamics and the Manna ei-state model. Results ob- 
tained here using FRG are shown for comparison. * It was 
computed using the scaling relation D = d + £ and the re- 
ported value of 
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